{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Paleopiezometry based on dynamically recrystallized grain size"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "module plot imported\nmodule averages imported\nmodule stereology imported\nmodule piezometers imported\nmodule template imported\n\n======================================================================================\nWelcome to GrainSizeTools script\n======================================================================================\nA free open-source cross-platform script to visualize and characterize grain size\npopulation and estimate differential stress via paleopizometers.\n\nVersion: v3.0.2 (2020-12-30)\nDocumentation: https://marcoalopez.github.io/GrainSizeTools/\n\nType get.functions_list() to get a list of the main methods\n\n"
     ]
    }
   ],
   "source": [
    "# Load the script first (change the path to GrainSizeTools_script.py accordingly!)\n",
    "%run C:/Users/marco/Documents/GitHub/GrainSizeTools/grain_size_tools/GrainSizeTools_script.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The script includes a function for estimating differential stress based on \"average\" recrystallized grain sizes named ``calc_diffstress()``. This function requires\n",
    "\n",
    "- defining the mineral phase and the piezometer relation to use,\n",
    "\n",
    "- entering the (apparent) grain size as the **equivalent circular diameter in microns**,\n",
    "- measured with a **specific type of \"average\" with no stereological correction**,\n",
    "- and set the type of stress, either uniaxial compression/extension or plane stress, for proper stress correction.\n",
    "\n",
    "For the first requirement, the GrainSizeTools script includes common mineral phases such as quartz, calcite, olivine and albite (more available soon) and a large list of piezometer relations (17 so far). Also, the script facilitates to write ad hoc piezometric relations.\n",
    "\n",
    "For the second requirement, the function will automatically convert the equivalent circular diameter to linear intercepts where applicable using de Hoff and Rhines (1968) correction. This is, **you don't have to worry about whether the piezometer was originally calibrated using linear intercepts**, always use the equivalent circular diameters in microns.\n",
    "\n",
    "The third requirement is key for a correct estimation of the differential stress since each paleopiezometer was calibrated for a specific average grain size (e.g. the arithmetic mean, median or RMS mean) and, hence, **only provides valid results if the same type of average is used**. Also, **you should not use any type of stereological correction for the estimation of the average grain size**, if the author(s) of the piezometer used any type of stereological correction, the average grain size will be automatically corrected by the function. \n",
    "\n",
    "The fourth requirement means that the user has to decide whether to correct or not the estimate of the differential stress for plane stress using the correction factor proposed by Paterson and Olgaard (2000). The rationale for this is that the experiments designed to calibrate piezometers are mainly performed in uniaxial compression while natural shear zones approximately behave as plane stress volumes.\n",
    "\n",
    "In the next subsection, we will show examples of how to obtain information about the different piezometers and define these parameters, but let's first load the script."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load the script first (change the path to GrainSizeTools_script.py accordingly!)\n",
    "%run C:/Users/marco/Documents/GitHub/GrainSizeTools/grain_size_tools/GrainSizeTools_script.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Get information on piezometric relations\n",
    "\n",
    "You can get information from the console on the different available piezometric relations  just by typing ``piezometers.*()``, where * is the mineral phase, either ``quartz``, ``calcite``, ``olivine``, or ``feldspar``. For example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "Available piezometers:\n'Cross'\n'Cross_hr'\n'Holyoke'\n'Holyoke_BLG'\n'Shimizu'\n'Stipp_Tullis'\n'Stipp_Tullis_BLG'\n'Twiss'\n"
     ]
    }
   ],
   "source": [
    "piezometers.quartz()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you want to get the details of a particular piezometric relationship you can do so as follows. Remember that the relationship between recrystallized grain size and differential stress is $\\sigma_d = B g^{-m}$ where $\\sigma_d$ and $g$ are the differential stress and the average grain size respectively."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": [
       "(550,\n",
       " 0.68,\n",
       " 'Ensure that you entered the apparent grain size as the arithmetic mean grain size',\n",
       " True,\n",
       " 1.5)"
      ]
     },
     "metadata": {},
     "execution_count": 4
    }
   ],
   "source": [
    "piezometers.quartz('Twiss')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note the five different outputs separated by commas which correspond with:\n",
    "- the constant *B* of the piezometric relation\n",
    "- the exponent *m* of the piezometric relation\n",
    "- A warning indicating the average to use with this piezometric relation\n",
    "- An indication of whether the piezometric relation was calibrated using linear intercepts (if ``False`` the piezometric relation was calibrated using equivalent circular diameters.\n",
    "- The stereological correction factor used (if applicable). If ``False``, no stereological correction applies.\n",
    "\n",
    "## Estimate differential stress using the ``calc_diffstress()`` function\n",
    "\n",
    "Let us first look at the documentation of the:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "text": [
      "\u001b[1;31mSignature:\u001b[0m \u001b[0mcalc_diffstress\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mgrain_size\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mphase\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mpiezometer\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mcorrection\u001b[0m\u001b[1;33m=\u001b[0m\u001b[1;32mFalse\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mDocstring:\u001b[0m\n",
      "Apply different piezometric relations to estimate the differential\n",
      "stress from average apparent grain sizes. The piezometric relation has\n",
      "the following general form:\n",
      "\n",
      "df = B * grain_size**-m\n",
      "\n",
      "where df is the differential stress in [MPa], B is an experimentally\n",
      "derived parameter in [MPa micron**m], grain_size is the aparent grain\n",
      "size in [microns], and m is an experimentally derived exponent.\n",
      "\n",
      "Parameters\n",
      "----------\n",
      "grain_size : positive scalar or array-like\n",
      "    the apparent grain size in microns\n",
      "\n",
      "phase : string {'quartz', 'olivine', 'calcite', or 'feldspar'}\n",
      "    the mineral phase\n",
      "\n",
      "piezometer : string\n",
      "    the piezometric relation\n",
      "\n",
      "correction : bool, default False\n",
      "    correct the stress values for plane stress (Paterson and Olgaard, 2000)\n",
      "\n",
      " References\n",
      "-----------\n",
      "Paterson and Olgaard (2000) https://doi.org/10.1016/S0191-8141(00)00042-0\n",
      "de Hoff and Rhines (1968) Quantitative Microscopy. Mcgraw-Hill. New York.\n",
      "\n",
      "Call functions\n",
      "--------------\n",
      "piezometers.quartz\n",
      "piezometers.olivine\n",
      "piezometers.calcite\n",
      "piezometers.albite\n",
      "\n",
      "Assumptions\n",
      "-----------\n",
      "- Independence of temperature (excepting Shimizu piezometer), total strain,\n",
      "flow stress, and water content.\n",
      "- Recrystallized grains are equidimensional or close to equidimensional when\n",
      "using a single section.\n",
      "- The piezometer relations requires entering the grain size as \"average\"\n",
      "apparent grain size values calculated using equivalent circular diameters\n",
      "(ECD) with no stereological correction. See documentation for more details.\n",
      "- When required, the grain size value will be converted from ECD to linear\n",
      "intercept (LI) using a correction factor based on de Hoff and Rhines (1968):\n",
      "LI = (correction factor / sqrt(4/pi)) * ECD\n",
      "- Stress estimates can be corrected from uniaxial compression (experiments)\n",
      "to plane strain (nature) multiplying the paleopiezometer by 2/sqrt(3)\n",
      "(Paterson and Olgaard, 2000)\n",
      "\n",
      "Returns\n",
      "-------\n",
      "The differential stress in MPa (a float)\n",
      "\u001b[1;31mFile:\u001b[0m      c:\\users\\marco\\documents\\github\\grainsizetools\\grain_size_tools\\grainsizetools_script.py\n",
      "\u001b[1;31mType:\u001b[0m      function\n"
     ],
     "name": "stdout"
    }
   ],
   "source": [
    "?calc_diffstress"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As indicated in the documentation, the ``calc_diffstress()`` requires three (obligatory) inputs: (1) the average grain size in microns, (2) the mineral phase, and (3) the piezometric relation to use. We provide a few examples below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "============================================================================\ndifferential stress = 83.65 MPa\n\nINFO:\nEnsure that you entered the apparent grain size as the arithmetic mean grain size\nECD was converted to linear intercepts using de Hoff and Rhines (1968) correction\n============================================================================\n"
     ]
    }
   ],
   "source": [
    "calc_diffstress(12, phase='quartz', piezometer='Twiss')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The function returns the differential stress (in MPa) plus some relevant information about the corrections made and the type of average expected. Most piezometric calibrations were calibrated using uniaxial compression deformation experiments while in nature most shear zones approximately behaves as plane stress. Due to this, it may be necessary to correct the differential stress value. The ``calc_diffstress()`` allows you to apply the correction proposed by Paterson and Olgaard (2000) for this as follows (note the slightly different value of differential stress):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "============================================================================\n",
      "differential stress = 96.59 MPa\n",
      "\n",
      "INFO:\n",
      "Ensure that you entered the apparent grain size as the arithmetic mean grain size\n",
      "ECD was converted to linear intercepts using de Hoff and Rhines (1968) correction\n",
      "============================================================================\n"
     ]
    }
   ],
   "source": [
    "calc_diffstress(12, phase='quartz', piezometer='Twiss', correction=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Some paleopiezometers require uncommon averages such as the root mean square or RMS, for example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": [
       "(669.0,\n",
       " 0.79,\n",
       " 'Ensure that you entered the apparent grain size as the root mean square (RMS)',\n",
       " False,\n",
       " False)"
      ]
     },
     "metadata": {},
     "execution_count": 8
    }
   ],
   "source": [
    "piezometers.quartz('Stipp_Tullis')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this case you should estimate the RMS as  \n",
    "$RMS = \\sqrt{\\dfrac{1}{n} (x_{1}^2 + x_{2}^2 + ... + x_{n}^2)}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "============================================================================\ndifferential stress = 36.79 MPa\n\nINFO:\nEnsure that you entered the apparent grain size as the root mean square (RMS)\n============================================================================\n"
     ]
    }
   ],
   "source": [
    "# Import the example dataset\n",
    "filepath = 'C:/Users/marco/Documents/GitHub/GrainSizeTools/grain_size_tools/DATA/data_set.txt'\n",
    "dataset = pd.read_csv(filepath, sep='\\t')\n",
    "dataset['diameters'] = 2 * np.sqrt(dataset['Area'] / np.pi)  # estimate ECD\n",
    "\n",
    "# estimate the root mean squared\n",
    "rms = np.sqrt(np.mean(dataset['diameters']**2))  # note that in Python the exponent operator is ** (as in Fortran) not ^ (as in Matlab)\n",
    "\n",
    "calc_diffstress(rms, phase='quartz', piezometer='Stipp_Tullis')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Estimation of the differential stress using arrays of values\n",
    "\n",
    "Alternatively, you can use (NumPy) arrays as input to estimate several differential stresses at once. In this case, the ``calc_diffstress()`` function will return a NumPy array, so it is generally more useful to store it in a variable as in the example below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "============================================================================\nINFO:\nEnsure that you entered the apparent grain size as the arithmetic mean in linear scale\nECD was converted to linear intercepts using de Hoff and Rhines (1968) correction\nDifferential stresses in MPa\n"
     ]
    },
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": [
       "array([167.41, 153.66, 162.16, 172.73, 162.83, 185.45])"
      ]
     },
     "metadata": {},
     "execution_count": 10
    }
   ],
   "source": [
    "ameans = np.array([12.23, 13.71, 12.76, 11.73, 12.69, 10.67])  # a set of average grain size values\n",
    "estimates = calc_diffstress(ameans, phase='olivine', piezometer='VanderWal_wet')\n",
    "estimates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If the set of estimated values belongs to the same structural element (e.g. different areas of the same mylonite or different rocks within the same shear zone), you may want to estimate the average differential stress from all the data. The GrainSizeTools script provides a method named ``conf_interval()`` for this."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "text": [
      "\u001b[1;31mSignature:\u001b[0m \u001b[0mconf_interval\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mdata\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mconfidence\u001b[0m\u001b[1;33m=\u001b[0m\u001b[1;36m0.95\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mDocstring:\u001b[0m\n",
      "Estimate the confidence interval using the t-distribution with n-1\n",
      "degrees of freedom t(n-1). This is the way to go when sample size is\n",
      "small (n < 30) and the standard deviation cannot be estimated accurately.\n",
      "For large datasets, the t-distribution approaches the normal distribution.\n",
      "\n",
      "Parameters\n",
      "----------\n",
      "data : array-like\n",
      "    the dataset\n",
      "\n",
      "confidence : float between 0 and 1, optional\n",
      "    the confidence interval, default = 0.95\n",
      "\n",
      "Assumptions\n",
      "-----------\n",
      "the data follows a normal or symmetric distrubution (when sample size\n",
      "is large)\n",
      "\n",
      "call_function(s)\n",
      "----------------\n",
      "Scipy's t.interval\n",
      "\n",
      "Returns\n",
      "-------\n",
      "the arithmetic mean, the error, and the limits of the confidence interval\n",
      "\u001b[1;31mFile:\u001b[0m      c:\\users\\marco\\documents\\github\\grainsizetools\\grain_size_tools\\grainsizetools_script.py\n",
      "\u001b[1;31mType:\u001b[0m      function\n"
     ],
     "name": "stdout"
    }
   ],
   "source": [
    "?conf_interval"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      " \n",
      "Mean = 167.37 ± 11.41\n",
      "Confidence set at 95.0 %\n",
      "Max / min = 178.79 / 155.96\n",
      "Coefficient of variation = ±6.8 %\n"
     ]
    },
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": [
       "(167.37333333333333, 11.412701448126, (155.96063188520733, 178.78603478145934))"
      ]
     },
     "metadata": {},
     "execution_count": 12
    }
   ],
   "source": [
    "conf_interval(estimates)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.9-final"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}